
// STL
#include <vector>     // std::vector<>

// jScience
#include "jstats.hpp" // mean() << jmm: temp >>


//========================================================================
//========================================================================
//
// NAME: double statistics::r(
//                            const std::vector<double> &x,
//                            const std::vector<double> &y
//                            )
//
// DESC: Computes the sample correlation coefficient.
//
// NOTES:
//
//     - See:
//
//     https://en.wikipedia.org/wiki/Pearson_product-moment_correlation_coefficient
//     http://mathworld.wolfram.com/CorrelationCoefficient.html
//
//========================================================================
//========================================================================
inline double statistics::r(
                            const std::vector<double> &x,
                            const std::vector<double> &y
                            )
{
    auto n = x.size();
    
    double x_mean = mean(x);
    double y_mean = mean(y);
    
    double num = 0.0;
    double x_denom2 = 0.0;
    double y_denom2 = 0.0;
    
    for(auto i = 0; i < n; ++i)
    {
        double xterm = (x[i] - x_mean);
        double yterm = (y[i] - y_mean);
        
        num      += (xterm*yterm);
        x_denom2 += (xterm*xterm);
        y_denom2 += (yterm*yterm);
    }
    
    return (num/(std::sqrt(x_denom2)*std::sqrt(y_denom2)));
}
